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ABSTRACT 

A direct numerical simulation of a turbulent channel flow with three passive scalars 
at different molecular Prandtl numbers is performed. Computed statistics including the 
turbulent Prandtl numbers are compared with existing experimental data. The computed 
fields are also examined to investigate the spatial structure of the scalar fields. The scalar 
fields are highly correlated with the streamwise velocity; the correlation coefficient between 
the temperature and the streamwise velocity is as high as 0.95 in the wall region. The 
joint probability distributions between the temperature and velocity fluctuations are also 
examined; they suggest that it might be possible to model the scalar fluxes in the wall 
region in a manner similar to the Reynolds stresses. 

INTRODUCTION 

The transport of heat and contaminants in turbulent flows is of great importance in 
many engineering applications. Heat-transfer problems in heat exchangers, gas turbines, 
and nuclear reactors, and pollution dispersal in urban atmospheres are a few examples. 
In general, a turbulent flow field causes fluctuations in a scalar field through turbulent 
convection whereas the fluctuating scalar field influences the velocity field through mean 
gradients and density changes. For small temperature differences or small concentration of 
contaminants, however, the turbulent velocity field drives the scalar field, and the influence 
of the latter on the former is rather weak and can be neglected. In these cases, a passive 
scalar field can be determined independently by solving the conservation equation of the 
passive scalar for a given turbulent velocity field. 

In solving the Reynolds-averaged conservation equation for a passive scalar, one en- 
counters the same turbulence closure problem as in solving the governing equations for the 
velocity field. In the simplest approach, one uses the effective eddy conductivity (or diffu- 
sivity) to relate the scalar fluxes to mean flow variables in analogy to the eddy viscosity in 
the momentum transport. Since the turbulent velocity field is assumed to be known, the 
eddy viscosity is known, and a knowledge of the turbulent Prandtl number (or Schmidt 
number for mass transfer), defined as a ratio of the eddy viscosity to the eddy conductivity 
(diffusivity), is sufficient to solve the conservation equation for the scalar field. A large 
number of experimental works on the turbulent Prandtl number for various fluids can be 
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found in the literature: see, Eckelman and Hanratty (1972), McEligot et al. (1976), Anto- 
nia (1980), Malhotra & Kang (1984), to name a few. Reynolds (1975) provides a survey 
on more than 30 models for predicting the turbulent Prandtl number. 

In the second-moment closure approach, one dispenses with the eddy-diffusivity model 
and directly solves the Reynolds-averaged transport equations for the scalar fluxes. This 
approach is equivalent to solving the transport equations for the Reynolds stresses in the 
momentum transport. In the second-moment closure, the transport equations for the 
scalar fluxes contain many terms that have to be modeled. But a justification for this 
approach is that approximations are made at more fundamental levels, and, hence, they 
will be applicable for more general situations. To some extent, this justification has been 
validated in the modeling of the momentum transport, and it is likely that this will hold for 
the transport of passive scalars. Interested readers should refer to Lumley & Khajeh-Nouri 
(1974) and Launder (1978) for further details. 

In this paper we study numerically the turbulent transport of passive scalars in a channel 
flow by directly solving the unsteady, three-dimensional Navier-Stokes equations and the 
equations for passive scalar fields simultaneously. The resulting solutions contain detailed 
information on the velocity and the passive scalar fields. Time-averaged statistics, as well 
as turbulence structures associated with the velocity field, are reported in detail by Kim 
et al. (1987). The objective of the present paper is to present results pertaining to the 
transport of passive scalars in turbulent flows. Other details such as budgets for the scalar 
fluxes, essential information for developing the second-moment closure models, will be 
reported elsewhere. 

Throughout the present paper, we use temperature as our passive scalar (heat transfer), 
but this can be replaced with a mass concentration (mass transfer). Likewise, for example, 
thermal diffusivity, heat flux, and Prandtl number can be replaced with, respectively, 
concentration diffusivity, mass flux, and Schmidt number. 

In the first part of the paper, turbulence statistics including the turbulent Prandtl 
number associated with fully developed scalar fields for three different molecular Prandtl 
numbers (0.1, 0.71 and 2) are presented. In the second part of the paper, we examine 
the correlations between the velocity and scalar fields. Air ( Pr = 0.71) was chosen as 
the medium for this purpose. The wall-layer structures identified by the temperature field 
are compared with the structures found in the velocity field. The role of the organized 
turbulence structures in scalar transport is discussed. 

NUMERICAL PROCEDURES 

The unsteady Navier-Stokes equations were solved numerically at a Reynolds number 
of 3300, based on the mean centerline velocity and the channel half-width (180 based 
on the wall shear velocity and the channel half-width), with about 2 x 10 6 grid points 
(128 x 129 x 128 in x, y, and z). A spectral method — Fourier series in the streamwise and 
spanwise directions and Chebyshev polynomial expansion in the normal direction — was 
used for spatial derivatives. The time advancement was made by a semi-implicit method: 
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the Crank-Nicolson scheme for viscous terms and the Adams-Bashforth method for the 
nonlinear terms. Further information on the numerical method can be found in Kim et al. 
(1987). No subgrid-scale model was used in the computation since the grid resolution was 
sufficiently fine to resolve all the essential turbulent scales. Grotzbach (1981) performed 
similar numerical simulations with subgrid scale models for high Reynolds number flows. 

Once the velocity field was advanced for each time step, the corresponding scalar fields 
were obtained by integrating the following conservation equations for scalars: 


dOi 

dt 


+ 


d 

dxj 


&iUj 


1 dHj 
Re Pri dxjdxj 


+ Qi 


(no summation over t) 


(1) 


where the $i (* = 1,2,3) represent three different scalars, Re and Pr t denote the Reynolds 
number and the molecular Prandtl numbers, and Qi represents a source term for each scalar 
field. All the variables are nondimensionalized by the wall shear velocity, u T , channel 
half-width, 6, and a reference temperature, 0 r . With these nondimensionalizations, the 
Reynolds number and the Prandtl number are defined as Re = UtS/v and Pr { = v/a t , 
where u and a, are kinematic viscosity and thermal diffusivity, respectively. 

The initial field for the velocity was taken from the simulation of Kim et al. (1987). This 
represents a fully developed turbulent channel flow. Two different initial and boundary 
conditions were used for the scalar fields. In case I, the initial and boundary conditions 
were given as follows: 


0,(x,j/,2,O) = 0.5(1 - y 2 ) 

6i{x, -1,«,0 = 0, 0 { (x, 1,2,*) =0 (2) 

with Qi = 2 1 (Re Pri). This represents a case in which the passive scalar is created 
internally and removed from both walls. In case II, the initial and boundary conditions 
were given as follows: 


0i(x,y,z, 0) = erfc( 
Oi(x,-l,z,t) = 1, 


1 + y x 
2y/RePri y 0 
0i(x,l,z,t) = 0 


( 3 ) 


where erfc denotes the complementary error function, and t/o is an arbitrary constant 
which determines the initial thickness of the thermal boundary layer at the lower wall. 
The source terms were turned off (Q,=0) for case II. This represents a case in which the 
passive scalars were introduced at the lower wall (y = -1), and removed from the upper 
wall (y = 1). In the present paper, all the results are from case I, unless otherwise specified. 
Details of the results obtained from case II will be reported elsewhere. 


STATISTICS OF PASSIVE SCALARS 
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The computations were carried out until the passive scalar fields reached statistically 
steady states. The initial velocity field was taken from a previous calculation in which the 
velocity field already had reached a statistically steady state. 

Profiles of the mean and rms fluctuations of the scalar fields from case I are shown in Fig. 
1. Since all the statistical quantities to be presented in this paper are symmetric about 
the centerline of the channel, only the statistics from one side of the channel are shown; it 
is implied that y is measured from that wall, although in some cases, y w is used explicitly 
to avoid confusion. The passive scalars are nondimensionalized by the wall surface scalar 
flux defined as 


o T = 


Qw 


1.90 

A/ r\ 

dy 


^ d $ ! 


a dy 

w 


pCpU T pCpU T 


U T 


where p and c p are the fluid density and the specific heat at constant pressure, respectively, 
and a = k/(pc p ). The surface heat flux 0 T is sometimes referred to as the friction tem- 
perature in analogy to the friction velocity, u T . Throughout the paper the superscript + 
indicates a nondimensional quantity scaled by the wall variables; for example, y + = yu T jv 
and Of = 0 x l0 Ti . 

In Fig. la, where the mean scalar profiles are shown in semilogarithmic coordinates, the 
formulas recommended by Kader (1981) are also included. Kader suggested this formula, 
which is a function of the molecular Prandtl number, as well as Reynolds number, for the 
mean temperature profile after examining several experimental results. As shown in the 
figure, the agreement between the computed results and Kader’s formula is remarkable. 
This is especially so since all the experimental data Kader considered were at much higher 
Reynolds numbers than the present one. 

The rms fluctuations for each scalar field are shown in Fig. lb, together with the exper- 
imental data from a slightly heated turbulent boundary layer obtained by Subramanian 
and Antonia (1981). The experimental data shown correspond to Pr = 0.71 and their 
lowest Reynolds number (Re m = 990 based on the momentum thickness and the channel 
half-width, compared to 280 for the present case). With the wall scaling, it is expected 
that the two results should agree in the wall region, and indeed the agreement is quite 
good. The locations of peak rms value move away from the wall as the molecular Prandtl 
number decreases. The thickness of thermal boundary layers defined as 


Ay = 



Oc, ~ dj 

°r. 


dy 


was 1.6, 2.4, and 2.8, for Pr = 0.1, 0.71, and 2, respectively, where 0 Ci is the mean 
temperature at the centerline. The corresponding Clauser thickness — a similar measure 
of the momentum layer — was 2.7. The variation in the thermal boundary-layer thicknesses 
indicates the breakdown of the “Reynolds analogy.” 

The normal heat fluxes, O'v ’ , and their correlation coefficients are shown in Fig. 2. With 
the present normalization, the normal heat flux represents a ratio of the normal heat flux 
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to the surface heat flux (6'v 1 /q w ). Again the symbols represent the data from Subramanian 
and Antonia. Note that, for the fully developed channel flow, —6'v' + (1 / RePr)6, y = 1 - y w 
when a steady state is reached, and the data are indeed very close to the straight line. The 
agreement with the experimental data is poor, but this is probably a result of the difference 
in Reynolds numbers and the difference between the two flows (channel vs boundary layer). 
For example, the correlation coefficient for the present case must be zero at the channel 
centerline, whereas the experimental data from the boundary layer can remain finite at the 
edge of the boundary layer. Figure 2b indicates that the correlation coefficient decreases 
with an increase in the molecular Prandtl number. 

Figure 3 shows the axial heat fluxes and their correlation coefficients. Again, the axial 
heat flux normalized in the present form represents a ratio of the axial heat flux to the 
surface heat flux. Unlike the normal heat fluxes, which should be < 1, the axial heat fluxes 
are > 1 for most regions except the lowest Prandtl number case. The agreement between 
the computed results and the experimental results is good away from the wall, but poor 
near the wall. The correlation coefficients between the temperature and the axial velocity 
are also much higher than those between the temperature and the normal velocity. In the 
wall region, the correlation coefficients are as large as 0.95 for the two high-Prandtl-number 
cases. It is also interesting to note that in the wall region, the correlation coefficients of 
6'u' for the higher molecular Prandtl number cases are larger than that for the lowest 
Prandtl number case, which is different from the trend of the correlation coefficients of 
6'v' . The agreement with the experimental data is satisfactory except near the centerline 
of the channel, where agreement is not expected. 

Profiles of eddy viscosity, diffusivity, and turbulent Prandtl number are shown in Fig. 
4. The low-Prandtl-number case is quite different from the other two, whereas there 
exists almost no difference between the cases of Pr — 0.71 and 2. The turbulent Prandtl 
number for the case of Pr = 0.1 remains > 1, and that for the case of Pr = 0.71 and 2 
has a maximum value of about 1.1 at the wall and decreases slightly away from the wall 
and remains at values < 1. From experimental observations, it is generally accepted that 
Pr T < 1 for Pr > 1, and Pr t > 1 for Pr < 1 (Reynolds, 1975), which is consistent with the 
present results. For all three cases, the turbulent Prandtl numbers have local maximums at 
y/6 «0.25 (y + ss45). It is interesting to note that the turbulent Prandtl numbers approach 
the same constant value (about 1.1) as the wall is approached, independent of the molecular 
Prandtl numbers. The same behavior was also observed for case II (not shown here). 
Deissler (1963) showed in his analysis of homogeneous turbulence with a uniform velocity 
gradient that the turbulent Prandtl number approached one at high velocity gradients 
independent of the molecular Prandtl number. 

STRUCTURE OF THE SCALAR FIELDS 

For case II, the time evolution of the scalar fields corresponding to temperature of 
air {Pr = 0.71) was also examined, since this made it possible to trace the turbulence 
structures originating from the wall region. Figure 5 shows contour plots of temperature 
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in an (x,y)-plane. The contours clearly show the existence of large-scale structures that 
bulge out from the wall region. Note that there exist sharp temperature gradients along 
the upstream edge of the large-scale structures. This structure was referred to as the 
temperature front by Chen & Blackwelder (1978) in their study of turbulence structures 
in a turbulent boundary layer. 

Contours of u' , 0 ', and O'v' in the (x,z)-plane at y + «5 are shown in Fig. 6. Contours 
of streamwise velocity fluctuations (Fig. 6a) display the well-known wall-layer streaky 
structures. Their mean spanwise spacing is about 100 in wall units. The contours of 
temperature fluctuations shown in Fig. 6b show nearly identical streaky structures. The 
regions of low and high temperature are elongated in the streamwise direction with a mean 
spanwise spacing of about 100 in wall units — the mean spanwise spacings were determined 
by examining the two-point correlations of the temperature fluctuations in the spanwise 
direction as in the case of the low-speed streaks (see Kim et al. 1987). It is clear from 
Fig. 6a and 6b that the temperature in the near-wall region is highly correlated with the 
streamwise velocity; the low-speed fluid is associated with low-temperature regions and 
the high-speed fluid with high temperature. In fact, the correlation coefficient at this y- 
location is about 0.95 (Fig. 3). Figures 6a and 6b are consistent with the experimental 
observation of Iritani et al. (1983), who performed simultaneous visualizations of velocity 
and temperature fields in a turbulent boundary layer, using hydrogen bubbles for the 
velocity field and surface-mounted, liquid-crystal sheets sensitive to the wall temperature 
for the scalar field. Figure 6c shows the contours of O'v'-, they display intermittent regions 
of large O'v' similar to the spatial distribution of u'v' . This implies that production of the 
scalar fluctuations (-O'v'dO / dy) also takes place intermittently just as that of the velocity 
fluctuations. In fact, contour plots of corresponding u'v' (not shown) are almost identical 
to those in Fig. 6c. 

Contours of the surface heat flux are shown in Fig. 7. By comparing this figure with Fig. 
6a we find that the regions of high surface-heat flux coincide with the regions of high-speed 
fluid. Regions covered by the low-speed streaks do not have high heat transfer at the wall. 
The rms fluctuation of the surface heat flux is 0.38 q w . The maximum surface heat flux 
in Fig. 7 is about 3 times the mean surface heat flux, and the minimum was about 40 % 
of the mean value, indicating that the distribution of the surface heat flux fluctuation is 
positively skewed. 

The high correlations between the streamwise velocity and temperature can be also seen 
in Fig. 8, where the joint probability distributions of ( u',v '), [u',0'), and [O' ,v') at ss 
5 are shown. Note that a straight line would indicate a perfect correlation, and Fig. 8 
indicates a strong positive correlation between u' and O ' , and a mild negative correlations 
between O' and v'. Because of the high correlation between u' and O', the joint probability 
distribution of O' and v' (Fig. 8b) is very similar to that of u' and v' (Fig. 8c). Examination 
of higher-order statistical correlations, such as skewness and flatness factors, also indicate 
a strong similarity between u' and O' in the wall region. This strong similarity in the 
behavior of the streamwise velocity and temperature can be observed also in the two-point 
correlations of the streamwise velocity and temperature fluctuations as shown in Fig. 9. At 
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this y-location, there exists practically no difference between R uu and Re t e , , even though 
the highest Prandtl number is 20 times the lowest one. The streaky structures of the 
streamwise velocity fluctuations, which can be identified with the long correlation length 
in the streamw'ise direction (order of Ax + « 1000) and the negative peaks in the spanwise 
correlation (at about Az + « 50), are so dominant in the wall region it appears that 
these structures are imposed onto the scalar fields almost independent of the molecular 
Prandtl number. This is rather surprising since the sublayer of the thermal boundary 
layer is a strong function of the Prandtl number (Fig. la). Obviously this Prandtl-number 
independence cannot be true in the limit Pr — * 0, and this is reflected somewhat in 
the two-point correlations for Pr = 0.1. This situation changes significantly away from 
the wall, where the two-point correlations of the temperature fluctuations depend on the 
molecular Prandtl number. The strong similarity between O'v' and u'v' has been reported 
by Perry and Hoffmann (1976), although their measurements were made at y/6 = 0.3 and 
the similarity between the two was not as strong as the present results. 

It appears that one can assume that the behavior of the passive scalar would be the same 
as that of the streamwise velocity fluctuation in the wall region. It should be interesting 
to examine whether one can model the fluxes of B'u' and O'v' in the same form as models 
for u' 2 and u'v'. 

SUMMARY 

Direct numerical simulation of a turbulent channel flow with three passive scalars was 
performed. No turbulence model was used in the computation, since the grid resolution was 
sufficiently fine to resolve all the essential turbulent scales. Computed statistics were com- 
pared with existing experimental data. The mean temperature profiles were in excellent 
agreements with the formula suggested by Kader (1981). Agreements in other quantities 
were satisfactory, and the general trend of the turbulent Prandtl number with respect to 
the molecular Prandtl number was consistent with the existing experimental observations. 
The computed scalar fields were also examined for structures in the scalar field. It was 
found that the temperature fields were highly correlated with the streamwise velocity. 
The streaky structures observed in the wall region were also observed in the temperature 
field with the same nondimensional mean spacing, almost independent of the molecular 
Prandtl numbers considered in the present paper. The correlation coefficient between the 
temperature and the streamwise velocity was as high as 0.95 in the wall region. The joint 
probability distributions between the temperature and velocity fluctuations and two-point 
correlations in the wall region suggest that it might be possible to model the scalar fluxes 
in a manner similar to that used for momentum fluxes. 
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Fig. 1 Profiles of a) mean and b) rms fluctuations of the scalar fields. 
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Fig. 2 Profiles of a) the normal heat flux, -0V and b) its correlation coefficient. 



Fig. 4 Profiles of a) eddy viscosity, b) eddy diffusivity and c) turbulent Prandtl number. 
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Fig. 5 Contours of a temperature field in an (x, y)-plane; mean flow direction is from left to right. 
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Fig. 6 Contours of a) streamwise velocity fluctuations, b) temperature, and c) 0V in the (x, *)-plane at 
y+ a 5. The solid lines represent negative quantities and the dashed lines represent positive quantities. 
The field shown here is a small portion of the computational domain, and the tick marks denote 50 i^/u r . 
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Fig. 7 Contours of the surface heat flux fluctuations. See Fig. 6 captions for legends. 
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Fig. 9 Two-point correlations of streamwise velocity 
and temperature fluctuations at y + ss = 5: a) stream- 
wise separation; b) spanwise separation. x/S — 1 corre- 
sponds tO XU T jv « 180. 


Fig. 8 The joint probability distributions of a) b) 0',t/',and 
c) u\v* at 5 for Pr^0.71. 
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